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1.0  INTRODUCTION 

The  computer  program  designated  ARO-1  (Ref.  1)  is  an  effective  tool  for  the  solution  of 
problems  definable  in  terms  of  the  three-dimensional,  unsteady  Euler  equations. 
Unprecedented  flexibility  and  applicability  of  the  code  to  a  broad  spectrum  of  aerodynamic 
flows  was  achieved  through  use  of  a  finite  volume  approach  in  Cartesian  coordinates.  The 
explicit  MacCormack  algorithm  expressed  in  fine-tuned,  CRAY-vectorized  coding  yielded  a 
very  fast  solver.  Unfortunately,  these  specific  attributes  which  make  the  code  so  useful  have 
delayed  correction  of  an  error  in  the  boundary  condition  routine. 

Proper  specification  of  boundary  conditions  is  a  difficult  problem  in  computational  fluid 
dynamics,  increasingly  so  with  added  dimensions  and  dependent  variables.  There  are 
actually  two  fundamental  problems  associated  with  boundary  conditions:  (1)  over-  or  under- 
specification  of  boundary  information  can  create  an  unstable  numerical  system,  and  (2) 
incorrect  specification  will  yield  invalid  results,  often  without  warning.  Cline  (Ref.  2), 
among  others,  has  demonstrated  that  the  theory  of  characteristics  can  be  used  to  deduce 
what  information  is  defined  in  the  interior  flow  and  thus  what  information  remains  to  be 
specified  at  the  boundary.  In  principle,  this  solves  the  first  problem.  The  second  problem 
requires  solution  through  experience  and  induction. 

Application  of  characteristics  theory  to  ARO-1  is  not  straightforward  because  of  the 
Cartesian  frame  of  reference.  This  report  presents  the  derivation  of  characteristic  relations 
in  generalized  coordinates  and  discusses  the  specialization  necessary  for  inclusion  in  ARO-1. 


2.0  COORDINATE  SYSTEMS  AND  NOTATION 


The  body-fitted  computational  grid  is  presumed  given  in  Cartesian  coordinates  (xJ,  x2, 
x3),  referenced  as  the  unbarred  system.  A  new  system,  termed  barred,  is  defined  locally  at  a 
boundary  point  as  sketched  in  Fig.  1.  The  basis  vectors  in  the  linear  barred  system  are 
constructed  in  the  following  manner.  Choose  ei  to  be  directed  from  a  boundary  point 
towards  the  associated  interior  point  such  that 


E|  •  n  =  1 

where  n*is  the  inward  unit  normal  to  the  boundary.  Select  e2  as 


ej  x  n 
|ej  x  n 


(1) 

(2) 
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Interior  Point 


Boundary 
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and  should  e,  x  n  =  0,  then  select  ^  as  a  unit  vector  in  an  arbitrary  direction  perpendicular 
ton.  Finally,  designate  £3  orthogonal  to  n and e2: 

e3  =  n  x  e2  (3) 

The  barred  system  is  linear  but  not  necessarily  Cartesian.  In  general,  e2  is  not  orthogonal  to 
ej.  The  x1  axis  is  directed  through  the  interior  point,  and  the  x2  and  x3  axes  lie  in  the  tangent 
plane. 

From  Eqs.  (1)  through  (3)  it  follows  that 


e,  .  e2  x  e3  =  1  (4) 

and  the  reciprocal  basis  vectors  are  given  by 

e1  =  e2  x  e3  =  n  (5) 

e2  =  e3  x  £1  =  e2  (6) 

e3  =  ej  x  O) 


It  should  be  noted  that  e1,  elf  e3,  and  e3  lie  in  the  same  plane;  that  e1,  e2,  and  e3  are 
orthonormal;  that  ej  and  e3  are  of  the  same  magnitude;  and  that 

ej  •  ek  =  6jk  (8) 


Since  both  the  barred  and  unbarred  systems  are  linear,  the  transformation  between  them 
is  expressible  as  a  matrix  of  constants, 


M*  =  -&L 

3  flxJ 


(9) 


consistent  with 


e3  =  Mj^efc 


(10) 


where  Mj<  is  the  kth  component  of  ej  when  expressed  in  terms  of  the  unbarred  basis  vectors. 
The  inverse  transformation 


m  = 


dxk  « 


is  obtained  consistent  with 


M?Wk  - 


(11) 

(12) 
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The  metric  tensors  are  given  by 

gjk  =  ej  •  ek  <13) 

gjk  =  ei  •  (14) 

where 

gja  gak  =  ^  (15) 

Note  that  the  definition  of  the  basic  vectors  implies  that  g11  =  L  Contravariant  and 
covariant  components  are  related  by 


iij  =  gjk  uk  (16) 

and  vectors  can  be  expressed  in  either  fashion;  for  example, 

q  =  ui  ij  =  Uj  ei  (17) 

The  magnitude  of  the  velocity  vector  is  given  by 

q  =  Vgjk  Uj  uk  =  Vgjk  uj  uk  =  VQj  fij  (18) 

For  subsequent  reference,  frequently  used  expressions  and  their  associated  inverses  are 
grouped  below: 


V'S 

II 

n  — , 

s 

4 

II 

&  a 

(19) 

ej  =  Mkek 

ej  =  Wfek 

(20) 

ek  =  Mfei 

^  =  Wkei 

(21) 

-  xk  =  M^xJ 

xk  =  Wf(xj  -  xi) 

(22) 

uk  =  Mkui 

uk  =  Wkui 

(23) 

Uj  =  Mkuk 

uj  =  Wkuk 

(24) 

The  “o”  subscripts  in  Eq.  (22)  refer  to  the  specific  boundary  point  which  is  the  origin  of  the 
barred  system.  Since  the  unbarred  system  is  Cartesian,  ej  =  eJ.  Scalars,  such  as  e  and  p, 
are  the  same  in  both  systems  and  are  not  barred. 
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Since  the  barred  system  is  linear,  the  metric  tensors  are  constant  throughout  space.  Thus, 
the  Christoffel  symbols  are  zero,  and  covariant  differentiation  is  simply  partial 
differentiation,  such  as 


uj.k  =  (25) 

In  spite  of  the  barred  system’s  having  been  completely  defined  here,  actual 
implementation  of  characteristic  boundary  conditions  does  not  require  computation  of  the 
second  and  third  basis  vectors  or  components. 


3.0  CHARACTERISTICS  OF  THE  EULER  EQUATIONS 
The  Euler  equations  can  be  written  as 


where 


Gt  +  F\k  =  0 


and 


Fk  = 


quk 


gu‘uk  +  gikp 
(E  +  p)uk 


and  i  and  k  assume  values  1,  2,  and  3.  The  equation  of  state  is  taken  to  be 

p  =  (y  -  1)  (E  -  qq2/2) 

These  equations  can  be  put  in  the  form 


where 


AQ,  +  Bk  Q-k  =  0 


(26) 

(27) 

(28) 

(29) 

(30) 

(31) 
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and 


1  0 

0 

A  = 

Uj  flj  Q 

0 

J-fl2  e5i 

1 

7  —  1_ 

uk 

0 

Bk  = 

u1  uk 

fij  euk  +  fij1 

cu1 

gik 

—  q2  uk 

2  4 

5jk(E  +  p)  + 

pUjUk 

i 

ESI 

7T 

(32) 


(33) 


with  i  and  j  assuming  values  1,  2,  and  3.  Equation  (30)  is  further  rearranged  to  yield 

A  Q,  +  B  Q-!  =  C  (34) 

where 

B  =  §'  (35) 

and 

C  =  -B^Qjp  -  B*Q-3  (36) 

Thus,  the  derivatives  with  respect  to  t  and  x1  are  kept  on  the  left-hand  side,  defining  a 
reference  plane  that  can  be  easily  treated  by  the  theory  of  characteristics;  everything  else  is 
transposed  to  the  right-hand  side. 

Let  X  and  T  be  an  eigenvalue  and  eigenvector  defined  by 

(X  A*  -  B*)  T  =  0  (37) 

Then,  multiplying  Eq.  (34)  by  T*  and  using  Eq.  (37),  one  obtains 

T*A(Qt  +  X  Q-,}  =  T*C  (38) 

Along  the  characteristic  direction  given  by 

dx1  =  Xdt 

dx2  =  0  (39) 

dx3  =  0 
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the  total  time  derivative  becomes 

+  («) 

and  Eq.  (38)  becomes 

T*A  =  T*C  (41) 

dt 

Equation  (41)  is  termed  the  compatibility  equation  and  gives  a  relation  among  the  total  time 
derivatives  in  a  characteristic  direction. 


Equation  (37)  has  five  eigenvalues,  only  three  being  distinct  since  one  of  them  is  triply 
repeated.  Appearing  in  the  solution  of  Eq.  (37)  is  the  speed  of  sound,  given  by 


a  =  yfyp/Q  (42) 

The  eigenvalues  and  eigenvectors  of  Eq.  (37)  are  as  follows: 

X,  =  u1  -  a  (43) 

T:-(>-l)(i-Q^+^47u',  -u,-  0/^,l)  <«> 

'  X2  =  u1  +  a  (45) 

T5=(Y-l)(iq^-7iTu'>-Si  +  Si'7iT,l)  W 

x3  =  u1  (47) 

t;  =  <7-0  (|-q!  -  -  a„l)  (46) 

X4  =  u>  (49) 

TJ-  j  (-u2,|n.o)  (50) 

X5  =  Qi  (51) 
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(52) 


4.0  ANALYTIC  SIMPLIFICATIONS 


The  left-hand  side  of  Eq.  (41)  involves  total  time  derivatives  in  a  characteristic  direction, 
whereas  the  right-hand  side  involves  spatial  partial  derivatives.  The  two  sides  will  be 
simplified  separately.  Two  of  the  factors  on  the  left-hand  side  can  be  multiplied  to  yield 


A 


dQ 

dt 


de 

dt 


dE 

dt 


(53) 


The  left-hand  sides  corresponding  to  the  five  eigenvectors  are 


T*  A 

dQ 

dt 

dp 

dt 

du1 
a  q  — 
dt 

(54) 

Tl  A 

dQ 

dt 

=  4p  + 

dt 

du1 

dt 

(55) 

t;  a 

dQ 

dt 

_  dp  _ 
dt 

a2  — 
dt 

4 

(56) 

TJ  A 

dQ 

dt 

du2 

dt 

(57) 

t;  a 

dQ 

dt 

duj 

dt 

(58) 

Now  the  right-hand  side,  Eq.  (36),  can  be  written  as 

C  =  -  («g  -  M“WpM£  Q-j  (59) 

Using  Eqs.  (11),  (23),  and  (31),  one  obtains 
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dxa 


W' 


3um 

dx° 


dp_ 

3x“ 


(60) 


Another  pair  of  products  can  be  reduced  using  Eqs.  (12),  (14),  (21),  (23),  (24),  and  (33)  to 
yield 


”  u* 

m?  q  0  ‘ 

Mf?  Bk  = 

Wj,un  ufl 

Wj,  Mj*  eu*  +  Mf  Wj,  quP  Wj, 

J_qV 

2 

M?(E  +  p)  +  Mj1  c  Ujj  7  up 

v-l 

(61) 


Let  Tj  be  the  components  of  T 

T  =  (ri)  (62) 

then,  from  Eqs.  (12),  (18),  (29),  (59),  (60),  and  (61) 

rc  —  («sr  -  m? w^){  tJ  ^  (eu^> 

+  ^i+l  W[j  (gUn  +  gnj3  p)  (63) 

+  ?s  ~  [(E  +  pV]] 

Define 

T*  =  lruri+l  Wi,,rs]  (64) 

then  Eq.  (63)  can  be  written 

T'C  =  -(ag  -  MfW^)  T*F^,a  (65) 


The  significant  simplification  here  is  that  the  local  coordinate  system  dependence  is  removed 
and  the  equations  are  expressed  in  a  global  Cartesian  system  consistent  with  ARO-1.  The 
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overbars  are  essentially  gone.  Listed  below  are  the  five  T  vectors  for  each  of  the 
corresponding  characteristic  directions. 


A 

Tf  = 

(7 

-  0  [^q2  +  ■ 

\  Q' 
y-  1 

.-Ui-Wiy.,1] 

(66) 

a 

n  = 

(> 

3-j-U1, 

-u,+  wy^.l] 

(67) 

A 

Tj  = 

(7 

a2 

-r 

“i.  l] 

(68) 

A 

t;  = 

_1_ 

e 

[-Uk.  gki.  0 

(69) 

II 

1 

e 

Mi 

[-Uk.  gki.  0 

(70) 

The  first  three  compatibility  equations  do  not  involve  u2  or  u3.  It  is  desirable  to  eliminate 
these  components  from  the  fourth  and  fifth  equations  as  well.  Combining  Eqs,  (24),  (41), 
(57),  (58),  (65),  (69),  and  (70),  the  last  two  compatibility  equations  can  be  written 

Ms =  -(53  -  MfW#  i  Ms  (~ Uk,  gki,  0)F?a  (71) 

with  s  =  2  and  3.  Define 


Ts-k  =  ^(-uk,  gki.  0) 

for  k  =  1,2,  and  3;  then  Eq.  (71)  can  be  written 


(72) 


From  differentiation  of  Eq.  (23)  one  obtains 


wi  duk  _  du1 
k  dt  dt 


(73) 


(74) 


Consider  the  system  of  three  equations,  Eq.  (74)  and  Eq.  (73)  with  s  =  2  and  3.  Note  that 
since  the  unbarred  system  is  Cartesian,  uk  =  uk,  and  the  system  can  be  written 
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r 

w| 

wi 

w] 

du1 

dt 

i _ _ 

Mj 

du2 

dt 

= 

m2(  Tg  c)  +  m|(t;c)  + 

M2(Tg  C) 

Mj 

Mj 

2 

MW 

l _ 

Lie. 

i _ 

m](tsc)  +  Mf(T$c)  + 

m3  (tj  c) 
— 

The  rows  of  the  coefficient  matrix  are  the  components  with  respect  to  the  unbarred 
coordinates  of  e1,  e2,  and  e3,  which  are  orthonormal.  Therefore,  the  coefficient  matrix  is 
orthogonal  and  its  inverse  is  its  transpose.  Thus,  one  can  directly  solve  Eq.  (75)  to  obtain 

T=w[|t  t;.„c  -  w;  E,  wj,(t;.„c)  to 

The  fourth  and  fifth  compatibility  equations  can  therefore  be  replaced  with  Eq.  (76)  with 
k  =  1,  2,  and  3,  which  does  not  involve  u2  or  u3  or  the  components  of  e2  or  e3.  This 
represents  a  significant  reduction  in  computational  burden.  Knowledge  of  the  tangent-plane 
variable  dependence  is  not  required  for  implementing  characteristic  boundary  conditions 
out  of  the  plane. 


5.0  FINITE-DIFFERENCE  APPROXIMATIONS 

All  of  the  right-hand  sides  of  the  compatibility  equations  involve  the  quantity 

C  =  -(df  -  Mf  W£)F0.a  (77) 

For  purposes  of  constructing  finite  differences  of  Eq.  (77),  three  pairs  of  points  are  chosen. 
Led  d|u  be  the  k’th  component  with  respect  to  the  unbarred  system  of  a  vector  from  the 
boundary  point  to  one  of  the  chosen  points,  a  indicating  the  pair  (1,  2,  or  3),  and  0 
indicating  which  point  in  the  pair  (1  or  2).  Let  d„  be  the  vector  from  the  first  point  of  the 
a-pair  to  the  second,  and  let  d&  be  its  components  with  respect  to  the  unbarred  system;  then 

d*  =  d*2  -  d*,  (78) 

The  three  pairs  of  points  are  chosen  so  that  the  vector  dj  is  in  the  same  direction  as  £|  and  the 
vectors  dj,  dj,  and  dj  are  not  coplanar.  Since  dj  and  e,  are  parallel,  let  c  be  their  constant  of 
proportionality;  then  by  Eq.  (9), 

df  =  c  (79) 
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Multiplying  by  yields 


c  =  df 

The  finite-difference  formula  is  derived  from 


(80) 


FJn0  =  F0  +  FJo,fd^  +  0(2)  (81) 

where  the  ct&  subscript  on  Fi  indicates  its  evaluation  at  that  point,  and  the  o  subscript 
indicates  its  value  at  the  boundary  point.  It  is  assumed  that  all  the  chosen  points  are  in  a 
neighborhood  of  the  boundary  point  where  the  second-order  terms  are  negligible.  Define 

FL  =  FL  -  Fii  ’  (82) 


Neglecting  second-order  terms,  Eqs.  (78),  (81),  and  (82)  yield 

Fj  =  FJ,,df 

a  otf  a 

Multiplying  by  D“  (defined  below),  one  obtains 

FL  Df  =  Htd'  Df 


(83) 


(84) 


Since  it  is  desired  to  obtain  an  approximation  to  Eq.  (77),  Df  is  defined  to  assure  that 

Dfd'  =  8 MjW]  (85) 

Selecting  the  three  pairs  of  points  such  that  the  da  vectors  are  not  coplanar  ensures  the 
existence  of  the  inverse,  df ,  such  that 


df  dk=  3k 

a  f  ua 


Multiplication  of  Eq.  (85)  by  dfand  using  Eqs.  (79)  and  (80)  yields 


Wj'Si 


By  using  Eqs.  (84)  and  (85),  an  approximation  to  Eq.  (77)  is  Finally  given  as 


C  =  -P.Df 


(86) 


(87) 


(88) 
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This  approximation  is  valid  in  a  small  region,  limited  to  the  neighborhood  where  the  second- 
order  terms  of  Eq.  (81)  are  negligible. 


The  procedure  for  computing  C  is  summarized  as  follows:  (1)  calculate  d*  and  their 
inverse  from  Eqs.  (78)  and  (86);  (2)  evaluate  Dj1  from  Eq.  (87)  in  which  W£  are  the 
components  of  n  expressed  in  the  unbarred  system;  (3)  compute  Fj,  with  Eq.  (82),  in  which 
Eq.  (28)  is  employed  in  the  unbarred  sense;  thence  (4)  C  is  given  by  Eq.  (88). 

i 

6.0  BOUNDARY  CONDITION  METHODOLOGY 

The  solution  is  known  at  time  t  =  tj,  and  it  is  desired  to  update  the  variables  on  the 
boundary  (x1  =  0)  at  t  =  t2-  Three  characteristic  lines  are  constructed  in  the  (t,  x1)  plane 
intersecting  the  boundary  as  sketched  in  Fig.  2.  The  analysis  can  be  divided  into  five  cases: 

a.  subsonic  inflow 

b.  supersonic  inflow 

c.  subsonic  outflow 

d.  supersonic  outflow 

e.  solid  wall 

For  boundary  condition  considerations,  the  terms  subsonic  and  supersonic  refer  to  the 
normal  component  of  the  velocity.  At  the  solid  wall  the  normal  component  is  zero,  and  only 
one  case  is  needed.  In  Fig.  2,  the  flow  is  known  at  the  points  marked  with  circles  since  these 
points  are  inside  the  computational  region  at  t  =  tt.  The  time  step 

At  =  t2  -  t!  (89) 

as  determined  by  the  Courant-Friedrich-Lewy  stability  criterion  in  ARO-1  is  such  that  the 
circled  points  are  on  the  x1  axis  between  the  boundary  point  and  the  associated  interior 
point.  Therefore,  the  values  at  points  designated  by  circles  can  be  obtained  by  interpolation 
rather  than  extrapolation.  The  flow  is  undefined  at  the  points  marked  with  triangles  since 
these  points  are  outside  the  computational  region  at  t  =  tj;  thus  the  compatibility  equations 
along  the  characteristic  lines  passing  through  the  triangles  cannot  be  used  and  must  be 
replaced  with  user-specified  boundary  conditions.  For  supersonic  inflow,  none  of  the 
characteristic  lines  can  be  used,  which  means  that  at  such  boundary  points  all  flow  variables 
must  be  specified.  For  supersonic  outflow,  all  of  the  characteristic  lines  are  used,  which 
means  that  nothing  can  be  specified;  the  flow  at  the  boundary  must  be  a  result  of  the  flow 
calculations.  The  remaining  three  cases  will  be  treated  separately. 
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a.  Subsonic  inflow 


c.  Subsonic  outflow 


b.  Supersonic  inflow 


e.  Solid  wall  (u1  =  0) 


Figure  2.  Types  of  boundary  conditions. 
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6.1  SUBSONIC  INFLOW 

For  subsonic  inflow  (Fig.  2a),  the  flow  transmits  information  to  the  boundary  by  way  of 
only  the  first  characteristic  line 


x»  =  (ui  -  a)t 


with  compatibility  equation 


dp 

dt 


(90) 


(91) 


This  equation  supplies  only  one  of  the  five  relations  needed  to  specify  the  boundary 
conditions.  For  subsonic  inflow,  the  stagnation  pressure  and  temperature  are  usually 
known,  and  the  assumption  of  isentropic  flow  is  generally  reasonable.  The  remaining  two 
conditions  are  a  constraint  on  the  direction  of  the  flow 

ui  =  (92) 

where  the  direction  cosines,  f j,  are  specified  at  the  boundary;  alternately,  a  zero  derivative 
could  be  specified  which  requires  the  flow  direction  at  the  boundary  to  be  the  same  as  at  the 
associated  interior  point. 


Reference  3  gives  isentropic  relations  for  pressure-velocity  dependence  as 


P  = 


and  sound  speed-velocity  as 


a  =  aT 


(93) 


(94) 


Rigorous  evaluation  of  these  expressions  would  require  expensive  iteration  at  each  boundary 
point  to  achieve  consistency  among  the  variables.  However,  linearization  in  time  avoids  that 
iteration  and  yields  consistency  in  the  following  fashion.  Differentiating  Eq.  (93)  and 
simplifying  with  Eqs.  (94)  and  (42),  one  obtains 

p  1  Aay 

Ap - Apr  -  oqAq  +  —  eq2 — —  (95) 

pr  '2  a| 
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Finite  differencing  Eq.  (91)  gives 

Ap  -  aeAO'  =  (T*  C)At  (96) 

Using  Eqs.  (23)  and  (92),  differentiation  yields 

Au1  =  f“Aq  +  qW^  At**  (97) 

Combining  Eqs.  (95),  (96),  and  (97)  results  in 

Q  [q  +  aW^l  Aq  =  —  ApT  +  QQ2-~-  -  a^q  A(°  -  (TIC)At  /gg) 

>-  J  pt  2  ay 

As  shown  in  Fig.  3,  the  point  designated  II  is  the  interior  point  that  was  used  to  define 
the  barred  coordinate  system,  evaluated  at  time  tj.  All  finite  differences  in  Eqs.  (95)  through 
(98)  are  taken  between  the  points  02,  the  boundary  at  time  ti,  and  Cl,  the  characteristic 
point  in  the  interior  at  time  t|.  Equation  (91)  is  valid  only  along  this  characteristic  direction. 


Figure  3.  Nomenclature  for  the  subsonic  inflow  case. 

In  summary,  the  procedure  for  updating  the  dependent  variables  Q  at  the  subsonic 
inflow  boundary  point  02  is  as  follows:  from  Eq.  (23), 

5J|  -  Wi  ug,  (99) 

where  W*,  are  the  components  of  IT  in  the  unbarred  Cartesian  coordinate  system. 

From  Fig.  3,  to  first-order  accuracy,  the  characteristic  slope  at  02  is  approximated  by 
that  at  01  to  yield 
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Xci  =  -  At(u^,  -  3^)  (100) 

From  Eq.  (22) 

x{,  =  Wi  (xft  -  xg,)  (101) 

and  interpolation  yields  the  dependent  variables  at  the  characteristic  point,  Cl,  as 


Qci  - 


The  local  Mach  number 


Qoi  +  -3p-(Qu 

*ii 


Kfc,  =  -SSL 


-  Qoi) 


aci 


is  used  to  evaluate  total  conditions  from  the  isentropic  relations  given  in  Ref.  3. 


aTci  =  aci  1  + 


[*  +  JL^-L(Mci)  ] 


PtCi  =  Pci  1  + 


7/(7  -  *) 


(102) 


(103) 


(104) 


(105) 


Since  total  conditions  are  prescribed  at  the  boundary  (possibly  as  a  function  of  time),  the 
differences  required  in  Eq.  (98)  are  given  by 

=  aTfl2  "  aTCj  (106) 

aPt  =  Ptoj  “  Ptc1  (107) 

Direction  cosines  are  computed  at  Cl  using 

«.  =  -HIl  (108) 

qci 

If  the  flow  direction  is  specified  at  the  boundary  (possibly  as  a  function  of  time),  then 

-  fg!  (109) 

or,  alternatively,  if  zero  derivatives  are  specified,  then 

=  0  (110) 
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and 


f02  “  ^Cl  (111) 

The  velocity  difference,  Aq,  along  the  characteristic  is  then  computed  from  Eq.  (98)  and  the 
velocity  on  the  boundary  is  given  by 


qo2  =  qci  +  Aq  (112) 

with  the  Cartesian  components  expressed  as 

uS2  =  qo2foa2  (113) 

State  variables,  p02  and  a02,  are  computed  using  the  isentropic  relations,  Eqs.  (93)  and  (94), 
of  Ref.  3  and  density  evaluated  from 


<?02 


7P02 

(a<)2)2 


(114) 


6.2  SUBSONIC  OUTFLOW 

For  subsonic  outflow  (Fig.  4),  the  interior  flow  transmits  information  via  the  first 
characteristic  line 


with  compatibility  equation 


x1  =  (u1  -  a)t 


dp  du1 

—  —  SLQ  -  =1  iC 

dt  dt 


and  along  the  coincident  third,  fourth,  and  fifth  characteristic  lines 

x1  =  u't 

with  compatibility  equations 


iE_(a)2^  =  T3C 

dt  dt 


_  T*/ 


and  from  Eq.  (76),  written  as 


■^T  =  +  ^  +  ”  Wi  -?t  W“^  +  aC 


(115) 

(116) 

(117) 

(118) 

(119) 
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with  k  =  1,  2,  and  3.  Although  Eq.  (119)  provides  three  equations,  considering  the 
dependence  of  u1  [see  Eqs.  (73)  and  (74)],  they  are  of  rank  two;  therefore,  one  further 
condition  remains  to  be  specified.  It  will  be  taken  as  a  specified  static  pressure,  p.  For 
practical  applications,  pressure  is  specified  at  one  isolated  boundary  point  and  zero  gradient 
in  the  outflow  direction  applied  at  all  other  boundary  points. 

Referring  to  Fig.  4,  the  solution  is  known  at  time  tj  at  the  points  01  and  II  and,  through 
interpolation,  at  the  characteristic  points  B1  and  Cl.  The  procedure  for  updating  the 
dependent  variables  Q  at  the  subsonic  outflow  boundary  point  02  is,  from  Eqs.  (22)  and 
(23): 

*ii  =  K  (xf,  -  *?,)  (120) 

^  =  w>8,  <121> 

where  are  the  components  of  n  in  the  unbarred  system.  From  Fig.  4,  with  first-order 
accuracy,  the  characteristic  slopes  are  evaluated  at  01  to  yield 

=  -Atu^  (122) 

x'ci  =  -At«  -  aoi)  (123) 


Figure  4.  Nomenclature  for  the  subsonic  outflow  case. 
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Interpolation  yields 


Qbi 

=  Qoi  +  (On 

Qoi) 

(124) 

xn 

Qci 

-  Qoi  +  %  - 

xn 

Qoi) 

(125) 

and  from  Eq.  (23) 

aBl  = 

(126) 

“Cl  =  WaUCI 

(127) 

A  finite-difference  representation 

of  Eq.  (116)  is  written 

as 

(p02  “  Pci) 

-  aciOci(iio2  “  ufcj) 

=  At[f*c]Ci 

(128) 

which  yields  u^2  where  p02  =  Pn  except  at  one  select  point  where  p02  is  user-specified. 

A  finite-difference  representation  of  Eq.  (118)  is  written  as 

(P02  -  Pbi)  -  (aBi)2(eo2  -  Obi)  =  Ai[t5c]b1  (129) 

which  yields  eoi-  Finally,  differencing  of  Eq.  (119)  can  be  constructed  as 

(u&  -  ug.)  -  Wi(u^  -  u{J  +  At  [tU*C  -  W>  B|  (130) 

which  is  solved  to  yield  uq2  which  then  completes  calculation  of  Q02  at  the  subsonic  outflow 
boundary. 

6.3  SOLID  WALL  BOUNDARY 

Information  from  the  flow  is  transmitted  to  a  solid  wall  boundary  (see  Fig.  2d)  via  the 
same  characteristic  lines  as  in  the  subsonic  outflow  case,  Eqs.  (115)  through  (119).  The 
additional  condition  needed  to  determine  the  flow  at  the  wall  is 

q  •  n  =  0  (131) 
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which  can  be  analyzed 


q  •  n  =  u“ea  •  e1  =  u“5^  =  0 


resulting  in 


Referring  to  Fig.  5,  the 


u1  =  0 

procedure  for  computing  Q  at  02  is  [from  Eq.  (22)1: 
*ii  =  -  xgj) 


032) 

(133) 


(134) 


Figure  5.  Nomenclature  for  solid  wall  boundary  conditions. 

From  Fig.  5,  assuming  the  characteristic  slope  at  02  is  approximately  that  at  01  (first-order 
accuracy), 

*ci  =  At  no,  (135) 

then  interpolation  yields 


and  from  Eq.  (23) 


Qci  =  Qoi  +  -^-(Qii  -  Qoj) 

*IL 


(136) 


*hi  -  W>ug, 


(137) 

t 
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A  finite-difference  representation  of  Eq.  (1 16)  is  written  as 

(P02  ~  Pci)  -  a<;ieci(uo2  -  Uci)  =  At[TJC]ci  (138) 

to  allow  calculation  of  P02  where  U02  =  0  from  Eq.  (133).  The  density  at  point  02  is  obtained 
from  a  finite-difference  representation  of  Eq.  (118)  in  the  form 

(P02  ~  Poi)  ~  (aoi)2  (eo2  ~  Coi)  =  At[TjC]oi  (139) 

and  the  Cartesian  velocity  components  evaluated  from 

3 

U&  =  uff,  +  At  [fj+aC  -  W'  g  W>T|+kc]01  (140) 

It  should  be  recognized  that  Eq.  (131)  can  be  relaxed  to  permit  a  specified  surface- 
normal  velocity  to  allow  imposition  of  an  auxiliary  wall  boundary  condition.  For  example, 
viscous  boundary-layer  effects  on  the  inviscid  Euler  solution  can  be  included  using  the 
concepts  developed  in  Ref.  4. 


7.0  CONCLUDING  REMARKS 

Characteristic  boundary  conditions  as  described  herein  have  been  implemented  in  the 
three-dimensional,  time-dependent  Euler  solver  known  as  ARO-1.  Experience  with  the 
modified  code  to  date  indicates  improved  accuracy  of  the  solutions,  particularly  for  internal 
duct  flow  computations  with  subsonic  inflow/outflow  conditions.  However,  the  accuracy 
improvement  relative  to  the  solid  wall  boundary  conditions  using  zero  normal  pressure 
gradient  is  significantly  less  than  expected.  Vorticity  is  generated  at  curved  surfaces  and 
manifested  as  total  pressure  losses  which  are  converted  streamwise.  Characteristic  boundary 
conditions  reduce  the  magnitude  of  these  losses  by  only  a  few  percent.  Computational  grid 
refinement  appears  to  be  necessary  to  minimize  the  total  pressure  loss,  but,  since  ARO-1  is 
an  explicit  code,  the  associated  increase  in  computer  time  is  generally  prohibitive. 

The  ARO-1  code  is  routinely  used  to  enhance  the  ground  test  capabilities  of  the  wind 
tunnels  and  engine  test  facilities  at  AEDC.  The  assumption  of  zero  normal  pressure  gradient 
as  a  wall  boundary  condition  has  proven  to  yield  sufficiently  accurate  solutions  for  most 
engineering  applications.  Whenever  increased  precision  is  deemed  necessary,  both  grid 
refinement  and  characteristic  wall  boundary  conditions  must  be  employed. 
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NOMENCLATURE 

A  See  Eq.  (32) 

a  Speed  of  sound,  Eq.  (42) 

B  See  Eq.  (35) 

Bk  See  Eq.  (33) 

C  See  Eq.  (77) 

C  See  Eq.  (36) 

c  Proportionality  constant  between  di  and  ei 

Dj  See  Eq.  (87) 

da  Vector  from  the  first  point  of  the  a-pair  to  the  second  point 

dk  kth  component  of  da;  see  Eq.  (78) 
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d„i3  kth  component  of  a  vector  from  the  boundary  point  to  the  /8-point  of  the 

a-pair 

d&  The  inverse  of  dk,  Eq.  (86) 

E  Specific  energy 

ej  ,  e2,  and  The  basis  vectors  for  the  barred  system,  Eqs.  (1)  through  (3) 
ej 

e1,  e2,  and  The  reciprocal  base  vectors  for  the  barred  system,  Eqs.  (5)  through  (7) 
e3 

Fk  See  Eq.  (28) 

Fo  The  value  of  Fk  at  the  boundary  point 

Fk  See  Eq.  (82) 

F&s  The  value  of  Fk  at  the  ,8-point  of  the  a-pair 

G  See  Eq.  (27) 

gjk  and  gik  Metric  tensors,  Eqs.  (13)  and  (14) 

f"  Direction  cosines  of  ql  Eq.  (92) 

M  Mach  number,  Eq.  (103) 

Mj1  Transformation  between  barred  and  unbarred  systems,  Eq.  (9) 

n  Unit  normal  at  the  given  boundary  point 

p  Static  pressure,  Eq.  (29) 

Q  SeeEq.  (31) 

q  Flow  speed,  Eq.  (18) 
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q  Velocity  vector 

•T  An  eigenvector 

Tk  The  kth  eigenvector 

T  See  Hq.  (64) 

Tk  See  Eqs.  (66)  through  (70),  and  Eq.  (72) 

t  Time 

Uj  and  ui  Covariant  and  contravariant  components  with  respect  to  the  barred  system 
of  q 

Transformation  between  the  barred  and  unbarred  systems,  Eq.  (11) 
xi  Barred  coordinates 

a  Dummy  repeated  index 

y  Ratio  of  specific  heats 

A  Finite  difference 

6*  Kronecker  delta 

X  An  eigenvalue 

Xt  The  kth  eigenvalue 

q  Density 

Tj  Components  of  T,  Eq.  (62) 

SUBSCRIPTS 

01,02,  Bl,  Points  shown  in  Figs.  3,  4,  and  5 
Cl,  and  II 
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T  Total  conditions 

o  Given  boundary  point 

SUPERSCRIPTS 

*  An  asterisk  indicates  the  transpose 

(  "  )  Overbar  indicates  association  with  the  local  boundary  coordinate  system 
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